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Abstract. Consider approximating a function / by an emulator f based on n observations of /. Let w be a point 
in the domain of /. The potential error of / at w is the supremum of \f(w) — g(w) \ among functions 
g that satisfy all constraints / is known to satisfy. The supremum over w of the potential error is 
the maximum potential error of /. Suppose / is in a known class of regular functions. Consider 
the set T of all functions in that class that agree with the n observations and are globally no less 
regular than those observations require / to be. We find a lower bound on the potential error of any 
emulator over / 6 T. This is the mini-minimax uncertainty. Its maximum over w lower-bounds the 
maximum potential error of any /. If that is large, every emulator based on these observations is 
potentially substantially incorrect. To guarantee higher accuracy would require stronger assumptions 
about the regularity of /. We lower-bound the number of observations required to ensure that some 
emulator based on those observations approximates all / 6 T to within e. For the Community 
Atmosphere Model, the maximum potential error based on a particular set of 1154 observations of 
/ is no smaller than the potential error based on a single observation of / at the centroid of the 
21-dimensional parameter space. We also find lower confidence bounds for quantiles of the potential 
error and the mean potential error over w in the domain of /. For the Community Atmosphere 
Model, these lower confidence bounds are an appreciable fraction of the maximum potential error. 
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1. Introduction. Emulators, also known as surrogate functions and metamodels, are im- 
portant tools for approximating functions that have been observed only partially. Common 
emulation methods include Kriging, Multivariate Adaptive Regression Splines (MARS), Pro- 
jection Pursuit Regression, Polynomial Chaos Expansions, as well as various Bayesian mod- 
eling techniques [1, 2, 3]. This paper studies the accuracy of emulators. 

A common use for emulators is to approximate computer models. A computer model 
can be considered a deterministic function / from a vector of inputs in the domain of / 
to an output. 1 Only a finite number of simulations can be conducted, though typically an 
intractable number of inputs are possible — for instance if any input parameter is a floating 
point number. By fitting an emulator to observations of /, it becomes possible to cheaply 
approximate the computer model. While many computer models are effectively black boxes, 
emulators typically have closed forms, amenable to analytic study. 
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J The function / might not be entirely deterministic, for instance, it could involve Monte Carlo simulations; 
moreover, in distributed parallel computations, numerical results can depend on the order in which subproblems 
happen to complete. These cases can be thought of as observing / + e, where e is stochastic noise. We do not 
address the complication of noise here; however, uncertainty in the observations of / only makes the problem of 
approximating / to a given level of accuracy more difficult. Because we focus on lower bounds on the difficulty 
of approximating / accurately, our results generally remain lower bounds when the observations of / are not 
only incomplete, but also imperfect. It is only necessary to add the step of finding a lower confidence bound 
on the regularity of / from the observations, rather than finding the lower bound directly. 
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Computer models known as HEB [4] may be particularly difficult to emulate: They de- 
pend on .High-dimensional inputs; they are Expensive to run; and they are effectively .Black 
boxes that are not amenable to closed-form, analytic study. Because the models are high- 
dimensional, it takes prohibitively many runs to explore the domain of / — the number of runs 
needed grows exponentially in the dimension. Because the models are expensive, actually per- 
forming so many runs is impractical or impossible. And because the models are black boxes, 
there are few (if any) constraints that can be exploited to extrapolate reliably from inputs 
actually tried to inputs not sampled. There are numerous practical HEB problems, such as: 

• Climate models: [5] (21-28-dimensional domains; 1154 simulations; Kriging and MARS) 

• Automobile crashes: [6] (15-dimensional domain; 55 simulations; polynomial response 
surfaces and artificial neural networks). 

• Chemical reactions: [7] (30-50-dimensional domain; boosted surrogate models) and 
[8] (46-dimensional domain; seconds per simulation). 

• Aircraft design: [9] (25-dimensional domain; 500 simulations; response surfaces and 
Kriging), [10] (22-dimensional domain; minutes per simulation; response surfaces and 
Kriging), and [11] (31-dimensional domain; 20 minutes to several days per simulation; 
Kriging). 

• Electric circuits: [12] (60-dimensional domain; 216 simulations; Kriging). 

But how accurate are emulators? For Bayesian emulators, it is common to use the posterior 
distribution to measure uncertainty [13]. Another common technique for "validation" is to 
measure the error of the emulator on observations that were not used to train the emulator [14]. 
Both techniques require conditions which usually either cannot be verified or are known to 
be false. Posterior distributions depend on the choice of prior and likelihood. However, 
input variables generally represent fixed parameters, e.g., physical constants, and thus are not 
random. 2 Even truly stochastic parameters seldom have known distributions. Validation using 
held-out observations is relevant only if the error of the emulator at the held-out observations is 
representative of the error everywhere, or at least at the unobserved inputs that might matter. 
If the observations were independent and identically distributed (HD), representative samples 
could be obtained through random sampling. But values of / are clearly not IID: nearby 
inputs produce similar outputs, at least when the computer model has some regularity. 3 

Of course, without any conditions on /, there is no basis for extrapolating from the values 
of / at observed inputs to the values of / at unobserved inputs. Even though / may be "a 
black box", it is often taken for granted that / is regular, in an unspecified sense. There are 
many ways regularity could be measured. If / is effectively a black box, there can be little 
scientific basis for measuring regularity in one way rather than another. 

We work with what is known in numerical analysis as the absolute condition number 
and in function approximation as the Lipschitz constant. If / has Lipschitz constant K, 
then the difference between the output at two input points is not more than K times the 

2 By "random," we mean stochastic. Bayesians use 'probability' to quantify uncertainties that are aleatoric 
(caused by randomness) and epistemic (caused by ignorance) . Using probability to quantify one's beliefs about 
a parameter does not make the parameter's value random. 

3 One might try to correct for the correlation between nearby inputs, a well-known problem (for instance, 
[15]). But adjusting for correlation requires assuming that the observations follow a known distribution — the 
same problematic assumption that a fully Bayesian approach requires. 
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distance between those two points (measuring distance in some pre-specified metric). Similar 
results could be derived for other measures of regularity, but Lipschitz bounds are particularly 
amenable to analysis. 

With perfect and complete knowledge of /, we could emulate it perfectly — by / itself. 
However, if the emulator / is constrained to be computable from the observations, without 
relying on any other information about /, accuracy cannot be guaranteed. If it were known 
that the Lipschitz constant of / is K, that, together with the observations, would make possible 
a guarantee of some level of accuracy. All else equal, the larger K is, the more difficult it is 
to guarantee that an approximation of / is accurate. 

The observations impose a lower bound on K (but no upper bound). We explore whether 
some approximation of /, computable from the observations alone, is guaranteed to be accurate 
throughout the domain of / — no matter what / is — provided / agrees with the observations 
and has a Lipschitz constant not greater than the observed lower bound on K. Viewed as a 
function of w in the domain of /, the minimax error of emulators over the set J- of functions 
that agree with the observations and have Lipschitz constant no greater than the lower bound 
is called the potential error. 

The potential error is the "mini-minimax uncertainty" the title refers to. The first "mini" 
refers to the regularity condition: since / is essentially guaranteed to be less regular than 
the n observations require, the potential error is a lower bound on the minimax uncertainty 
for functions that are as regular as /. The second "mini" refers to emulators: this is the 
uncertainty for the best emulator. The "max" is over functions that agree with / at the 
observations and satisfy the optimistic regularity condition. That is, the potential error is the 
smallest that the maximum uncertainty could be, for the best emulator, for functions no less 
regular than / is required to be by the observations. The maximum potential error over the 
domain of / is the maximum potential error. 

If K were known, finding the potential error would be a standard problem in information- 
based complexity [16, 17, 18]. However, K is unknown because / is only partially observed. We 
derive bounds on the potential error using a lower bound for K computed from the observed 
variation in /. 

First, we derive a lower bound on the number of additional observations that could be nec- 
essary to learn / (section 3). We apply the bound to climate simulations from the Community 
Atmosphere Model, and find that the number of computer simulations needed to approximate 
/ to within any useful accuracy could be astronomical. 

Next, we give two lower bounds on the maximum potential error for approximating / from 
a fixed set of observations (section 4). The first bound is empirical. If this bound is large, 
then every emulator based on these observations is potentially substantially incorrect. The 
second bound, which is tighter, is a fraction of the unknown Lipschitz constant. That bound 
enables us to give conditions under which the constant emulator based on just one observation 
at the centroid has smaller maximum potential error than any emulator trained on the actual 
observations. When these conditions hold — and they hold for the climate simulations — every 
emulator is potentially substantially incorrect, even though there might be many observations, 
made at great expense. 

We then extend the results for the maximum pointwise error over the domain of / to 
quantiles of the error across the domain, and to the mean of the error across the domain. 
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Finally, we consider additional information about / that would reduce the potential error, 
and the availability of such information in applications. We also allude to measures that 
might be more relevant than potential error in some scientific or engineering contexts. 

2. Notation and problem formulation. The function / to be emulated is a fixed unknown 
real- valued function on [0, l] p , the p-dimensional unit cube. The space of real- valued contin- 
uous functions on [0, 1} P is C[0, 1} P . The Roman letters i, j, n, N, p, and q denote integers. 
Lowercase Greek letters denote real scalars, with the exception of [i, which denotes Lebesgue 
measure. Uppercase Roman letters such as X and D denote subsets of [0, l] p ; X is a fixed 
finite subset of [0, l] p . Lowercase Roman letters from the end of the alphabet, such as v , 
w, x, y, and z, denote points in [0, l] p . The lowercase Roman letters e, /, g, and h denote 
real- valued functions on (subsets of) [0,l] p . The domain of the function g is dom(<7). The 
restriction of g to D C dom(g) is denoted g\o- The observations from which / is to be emu- 
lated are f\x- we observe / on X. (We assume that the observations are noise- free; section 6 
discusses adapting the method to noisy observations.) The function / denotes a function on 
[0, l] p that may be selected using the observations f\x, but no other information about /; it is 
an emulator of /. Let ||/i||oo = sup^^m^) the infinity-norm of h. This paper concerns 

how large ||/ — /||oo could be, for the best possible choice of /: the minimax pointwise error 
among emulators. 

Let d be a metric on dom(g). The (best) Lipschitz constant for g is 

K(g) = sup \ ^^ V \ — ^ : v, w € dom(g) and v ^ wX . (2.1) 
I d(v,w) J 

If / £ C[0, If, then K(f) = oo. Define 

F K (g) = {hE C[Q, If : K(h) < k and h\ dom{g) = g}. 

Then J : o{f\x) is the space of functions in C[0, l] p that fit the observations: they interpolate 
/ on X. 

Definition. The potential error of f £ C[0, 1} P over the set of functions J- is 

= sup 1 1 f{w) - g(w)\ : g G J"| . 
Definition. The maximum potential error of f € C[0, l] p over the set of functions J- is 
£(f,T)= sup £(w;f,7) = hf-g\\ 00 :g€7\. 

wG[0,1]p 1 J 

Maximum potential error is an instance of worst-case error, as defined in information- 
based complexity literature [16, 17, 18]. The uncertainty of the emulator / can be quantified 
by its maximum potential error over the set J'ooiflx) of continuous functions that agree with 
the observations, £ (/, ^*oo(/|x)- 

This measure of uncertainty presumes that / G C[0, l] p . If / ^ C[0, l] p , f could differ 
from / by even more. Our main results concern lower bounds on the uncertainty of the best 
possible emulator of /, under optimistic assumptions about the regularity of /. Why make 
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Figure 2.1. The dotted line is tangent to f where f attains its Lipschitz constant. Thus it has slope K . The 
dashed line is the steepest line that intersects any pair of observations. Thus it has slope K . Clearly K < K . 



assumptions about the regularity of /? The set X is not dense in [0, l] p , so for any c > 0, 
there exists some function g 6 T 00 {f\x ) satisfying ||/ — g^oo > c - Thus £ (f, Fooiflx)^ = oo: 
the maximum potential error is infinite unless / has more regularity than mere continuity. 

Let K = K(f) and K = K(f\ x )- Because X C [0, K < K. Figure 2.1 illustrates this 
inequality. In this and subsequent figures, p = 1 and the bold black dots represent f\x, the 
observations of / at x G X. 

Henceforth, to simplify notation, we set 



X 



and 



Define the radius of T C C[0, 1} P to be 



:g,heT}. 



By the triangle inequality, for any /, 

£ K (f) > r{T K ). 

Equality holds for the emulator that "splits the difference": 



inf g(w) + sup g(w) 



(2.2) 



[16, 17]. That is, for all emulators / that agree with / on X, 



£«(/) > £M = £* K 
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Figure 2.2. For both panels, K = and the optimal interpolant f* is constant. In the left panel k — K . 
In the right panel k < K. If k, > K then e«T < / < ej, and, equivalently, f € T K . 




Figure 2.3. The length of the black error bars are twice the maximum potential error over J- K . The 
succession of panels shows that as the slope between observations approaches k, e*(w) approaches for points 
w between observations, and the maximum potential error over J- K decreases. 

the emulator /* is a minimax (over / € J- K ) for infinity-norm error. The function /* depends 
implicitly on the observational design X through T K = T K {f\x)- Because K < K, 

J 7 ^ C Tk 

and hence 

E* K > S* k . (2.3) 

That is, the minimax maximum potential error over the set of functions that agree with / 
on X and have Lipschitz constant K is at least as large as the minimax maximum potential 
error over the set of functions that agree with / on X and have Lipschitz constant equal to 
the empirical Lipschitz constant K — the least regularity consistent with the observations f\x- 
Define 

e+(w) = et x (w) = mm [f(x) + nd(x, w)] , 



Uncertainty Quantification for Emulators 



7 



and 



The last of these, e^(w), measures minimax error at w: Consider the set of functions g that 
agree with the observations f\x and have and have Lipschitz constant no larger than k. If 
K < K, for all we know / could be any element of that set. Consider all possible emulators 
/ that also agree with the observations f\x- The smallest (across emulators /) maximum 
(across functions g) error at the point w E [0, l] p is e*(io). Figures 2.2 and 2.3 illustrate these 
definitions. 

3. Bounds on the number of observations needed to approximate / well . To determine 
how many observations are required to approximate /, we must first decide how much error 
e > to tolerate in the approximation. 

Definition. If /\a — g\A < e> then / e- approximates g on A. If A = dom(g), then / 

oo 

e- approximates g. 

Definition. If J 7 is a non-empty class of functions with common domain D, then / e- 
approximates T on A C D if Mg £ J 7 , / e- approximates g on A. If A = D, then / 
e- approximates T . 

The emulator / e-approximates J- if and only if the maximum potential error of / on J- 
does not exceed e. 

If ultimately we seek to e-approximate what is a reasonable value for e? Since K is 
the observed variation of / on X, a useful value of e would typically be much smaller than 
K. (Otherwise, we might just as well take / to be a constant; c.f. section 4.) For fixed / and 
e, the number of observations needed to ensure that /* e-approximates Tk depends on the 
experimental design X. The following definitions help to bound the number of observations 
needed, without restricting our analysis to a particular experimental design. 
Definition. For fixed e > 0, and Y C dom(/), Y is e-adequate for / on A if /j^ e-approximates 
FnifW) on A. If A = dom(/), then Y is e-adequate for /. 

Let B(x,5) denote the open ball in W centered at x with radius 5. Let 

Nj = min{#y : Y is e-adequate for /}, 
where #Y is the cardinality of Y. 

Definition. The minimum potential computational burden is 

M = max{N g : g G J~k}- 

This is also a minimax quantity: over all experimental designs Y, M is the smallest number 
of observations of / needed to guarantee that the maximum error of the best emulator based 
on those observations is not larger than e. Minimum potential computational burden is an 
example of minimum worst-case cost, as defined in information-based complexity literature 
[17, 18]. 
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3.1. An upper bound on Nf. For each x € X, f K e-approximates TkHIk) on ( a t least) 
B(x,e/K). Thus, f K e-approximates Tk on IJ^gx B(x,e/K). Hence, the cardinality of any 
Y C [0, If for which 



V = {B (x, -1) : x € y} D [0,1] J 



is an upper bound on iVy. Because if is unknown, however, j^V has limited utility as an 
upper bound on Nf in practice. 

To generalize this formulation, let K + be a known upper bound on K. When d is such 
that [0, l] p may be covered by N + balls of radius e/K + , then Nf < N + . While, in general, 
finding an optimal covering set is NP-hard, it is often straightforward to find good upper 
bounds, or even exact solutions, for some metrics. For example, in i^, [0, l] p can be covered 



by 



K+ 
2e 



balls of radius e/K + . 



3.2. A lower bound on M. Depending on / and X, it can happen that e-approximates 
Tk on regions of the domain not contained in U x ^xB(x,e/K). To see this, consider p = 1, 
f(x) = x, and X = {0, 1}. Then K = K = 1. In this case, the observations f\x determine 
/ exactly: the only function in Tk is /• In this example, for a function g to agree with the 
observations requires it to attain the Lipschitz constant K everywhere. A function cannot 
agree with the observations and "run away" from /. 

More generally, if / varies on X, then for a function g to agree with / at the observations 
requires g to vary too. That required variation "spends" some of g's Lipschitz constant, 
preventing g from running as far away from / as it could if fx were constant. We now 
quantify this intuition to find lower bounds for M. 

Define 7 = arg min 7eK X^eX \f( x ) — 7l P - Computing 7 is straightforward because the 
objective function is univariate and convex. 4 Let X + = {x € X : f(x) > 7} and let X~ = 
{x E X : f{x) < 7}. Let 



xex+ 



and 



xex- 



Then Q + n Q~ 



Alternatively, we could set 7 = X^teJt f( x )i wnere #X is the size of X. The resulting lower bound 
may not be as tight. 

5 Fix x + £ X + and x~ £ X' . Then | f{x + ) - f(x~)\ /d(x + , x~) < K. Equivalently, d(x + ,x~) > 

\f{x + ) - f(x~)\ /K. Let B + = B (x + , [f{x) -y]K\ and B~ =b(x~, [y - f(x)] /k). Let a be the sum of the 

radii of B+ and B~ . Then a = (/(a; + ) - 7) /K + (7-/(2;")) /if = - /(as - )) /A', and a < d(ai + ,ar). 

Therefore, -B + n B~ — 0. Because our selection of x + G X + and x~ £ X~ was arbitrary, Q + n Q~ = 0. 
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Figure 3.1. The function f (shown in the left panel) is comprised of segments of z\, and the constant 
function 7 (all shown in the right panel). It is constant over roughly half of the domain. No function between 
e~~ and e"t (inclusive) is constant over a larger fraction of the domain. 



Define 



/:[0,lf 



e - (w), w £Q+ 
w 1 y ^ e^t (w), w G Q- 



1, 



otherwise. 



Figure 3.1 illustrates this definition. If we know f\x, we know /. By construction, / £ C 
Tk- 

Let Q = [0, l] p \ (Q + U Q-). Let /i be Lebesgue measure. By the union bound, because 
M([0,1PO = 1, 

KQ) > 1 " E » { B (x,\f(x)-7\/K)). 

p/2 

Let C2 = r(p/ 2+1 ) an d Coo = 2 P , where T is the gamma function. Then, for q £ {2, 00}, 



m(q)>i-c ? E (1/0*0 -7l/£) • 

If there is some x £ X for which the singleton set {x} is e-adequate for / on A C Q, then 
^(j4) <fi(B(0,e/K)). Hence, because / € J 7 ^-, 



M > 



> 



fj,(B(0,e/K)) 



(3.1) 



This lower bound can be extremely large when p is not small. 
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3.3. Application to climate modeling. Under the Uncertainty Quantification Strategic 
Initiative Laboratory-Directed Research and Development Project 10-SI-013, the Lawrence 
Livermore National Laboratory (LLNL) Institutional Science and Technology Office pro- 
vided us with results from 1154 climate simulations using the Community Atmosphere Model 
(CAM). Each of the p = 21 parameters was scaled so that the interval [0, 1] contained all val- 
ues considered physically reasonable. The output of interest was a scalar, the simulated global 
average upwelling longwave flux (FLUT) approximately 50 years in the future. That is, the 
simulator could be thought of as a function / that maps [0, l] p — > M. Observing / by running 
the simulator for a particular input was computationally expensive; each run took several days 
on a supercomputer. LLNL used several approaches to choose the points X C [0, l] p at which 
to run simulations, including Latin hypercube, one-at-a-time, and random-walk multiple-one- 
at-a-time [5]. The 1154 simulations include all points selected by any of those approaches. 

For these observations, 7 = 232.77; K = 14.20 for q = 2. By (3.1), 

> e- 21 x 10 26 . 

For example, if e is 1% of K, then M > 10 43 . Even if e is 50% of K, M > 10 8 . For q = 00, 
K = 34.68; in that case 

> e~ 21 x 10 25 . 

These lower bounds on the minimum potential computational burden are extreme for a 
wide range of values of e: there are functions that fit the observations and are no less regular 
than the observations require, but that cannot be approximated with useful precision from 
any tractable number of observations. The function /, which is simple to construct, attains 
these lower bounds on potential computational burden. 

4. Bounds on the maximum potential error for a fixed experimental design. The pre- 
vious section gave lower bounds on the potential computational burden (cardinality of X) 
required to attain a desired maximum potential error e. This section gives two lower bounds 
on the maximum potential error £^ for a fixed experimental design X: an absolute bound and 
a bound expressed as a fraction of K. The bound as a fraction of K can yield a strong negative 
result: when a statistic — calculable from the observations — exceeds a calculable threshold, the 
maximum potential error is no less than the maximum potential error from a single obser- 
vation at the centroid of the domain. That is, observing / for every point in X was (from 
the perspective of maximum potential error) wasteful: a single observation would have been 
better. In the CAM climate model example, maximum potential error is quite large, both in 
magnitude and fraction of K. 

4.1. Lower bounds. For g G C[0, l] p , define 

inf0 = M{g(w) : w G [0, l] p } 

and 

sup g = sup {g(w) : w E [0, l] p } . 



M > 



-21 



1.57 x 10 
0.0038 



si 



6.81 x 10 



21 



M > 



-21 



2.19 x 10 



6.81 x 10 



25 
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Proposition 1. supe* = £*. 

The proof of this proposition and others is in A. 
Because C Tk, for any emulator /, 

£ K {f) > r{F k ). (4.1) 

Hence, by proposition 1, we have established we have established 
Corollory 2. £x{f) > supe*-,. 

Corollary 2 is one of our principal results: it shows that supe*^,, a statistic calculable solely 
from the observations f\x, is a lower bound on the maximum potential error for any emulator 
/ based on the observations f\x- 

Theorem 3 gives a stronger lower bound in terms of the unknown value of A. 
Theorem 3. For any A, if supe^. > AA', then £x{f) > A A. 

4.2. Maximum potential error for an emulator based on one observation. In sections 4.2 
and 4.3 we work in l^: d(v,w) = \\v — w\\oo, which makes the calculations simpler and gives 
a particularly strong result. 

Let z = (1/2, ... , 1/2), the centroid of [0, l] p . Let g € J z oo(f\{ z }) be the constant function 
g{w) = f(z), Vm S [0, l] p . The distance from z to any point on the boundary of [0, 1} P is 
1/2, so 

That is, the maximum potential error of the emulator that is constant throughout [0, l] p and 
equal to the value of / at the centroid of the cube is A/2. Let W C [0, l] p be finite and c € M. 
Suppose f\w = c Let h £ J~oo(f\w)- By examining the corners of the domain, it follows that 
if \W\ < 2 p , 

£ K (KT K (f\w)) > y- 

That is, if / is constant on W, any emulator based on fewer than 2 P observations of / will 
have at least A/2 maximum potential error. Making 2 P observations of / is intractable for 
CAM and many other applications. If / is nearly constant the situation may still be hopeless. 

How do we know whether f\x is too close to constant to benefit from observing it more 
than once, but fewer than 2 P times? 
Corollory 4. If sup e*^ > A/2, then 

£kU) = £kU^kU\x)) > § > e K (g,? K U\{z}))- 

That is, if supe*^ > A/2, no emulator based on observing f\x has smaller maximum 
potential error than the constant emulator based on a single observation — / is too nearly 
constant. Corollary 4 follows directly from theorem 3, taking A = A/2. 

4.3. Application to climate modeling. We now return to the Community Atmosphere 
Model. Is the maximum potential error of the best emulator based on observing / at the 1154 
points in X lower than the maximum potential error of the constant emulator based on one 
observation of / at the centroid of [0, l] p ? We cannot simply compute these two maximum 
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potential errors, because K is unknown. But corollary 4 applies if we can determine whether 
supe^. > K/2. Recall that K = 34.68 in for this dataset. Unfortunately, determining 
sup is difficult. In i^, if f\x is constant, it amounts to finding a maximal empty hypercube, 
a problem recently shown to be NP-hard in p [19]. It is generally no easier if / varies on X. 
Fortunately, it suffices to bound supe*-.. Since we are working in l^, we can bound supe^, 
above and below by considering just the corners of [0, l] p . 

4.3.1. Computing an upper bound from non-adjacent corners in l^. Recall that through- 
out sections 4.2 and 4.3 d(v,w) = \\v — w\\oo. 

Lemma 5. Let = (0, ... ,0) and 1 = (1, . . . , 1). For v,w G [0, If, 

d(v, w) < d{v) = max (d(v, 0),d(v, 1)) . 

Proposition 6. 

supe*-. < - jmin f(x) + Kd(x) — max f(x) — Kd[x) \ . 

2 I x&X L J xgX L J J 

Using this proposition, we calculate supe*^ < 20.95 for the CAM dataset. 

4.3.2. Computing lower bounds from corners in l^. Clearly 

supe^ > maxje^(w) : Vw G {0, 1} P | . 

Perhaps surprisingly, this lower bound is essentially sharp for the CAM dataset. The domain, 
[0, l] p , contains 2 P corners {ri}f . Divide [0, l] p into 2 P hypercubes {-Rj}|f :1 with edge-length 
1/2, disjoint interiors, each containing a different corner of [0, l] p (e.g., one such hypercube 
is [0, l/2] p ). Then the R{ are disjoint ^-balls of radius 1/2. Because X contains only 1154 
points, the vast majority of the i?j do not contain any element of X. Because e*^ tends 
to increase with distance from points in X, these unoccupied hypercubes are good regions 
to explore to find points maximizing e^.. Within an unoccupied hypercube Ri, no point is 
farther in from any point in X than the corner fj. So, the corners {ri}f are good places 
to observe e*> for the purposes of establishing a tight lower bound on supe*>. 

For the CAM dataset, one corner r 7 - attains e*- (r,) = 20.95. So, e*~ attains the upper 
bound established in the previous section, and we conclude that supe^. = 20.95. 

4.3.3. Implications for the Community Atmosphere Model. Because supe*-, = 20.95 > 

17.34 = K/2, theorem 3 says that Sxif) ^ K/2 for any interpolation /. In other words, by 
the discussion in section 4.2, our maximum potential error would have been no greater had 
we just observed / once, at z, and predicted f{w) = f(z) for all w G [0, l] p . 

In some sense, this result is not surprising: if we had fixed K but replaced / with a constant 
function, and j^X < 2 P , then supe*-, > K/2, with equality holding if and only if z E X. By 

repeating the bounding procedures from the previous two sections with K/2 = 17.34 fixed but 
/ replaced with constant function c, we find supe*^,^ = 26.95. The increase in maximum 
potential error from 20.95 to 26.95 that results from replacing / with a constant shows that 
the observed variation in / reduces the maximum potential error considerably — although the 
maximum potential error remains quite large. 
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95% 


lower confidence bound 




norm 


lower quartile 


median upper quartile 


average 


Euclidean 


1.454 


1.596 1.731 


1.595 


supremum 


0.649 


0.717 0.782 


0.715 



Table 5.1 

Error of the minimax emulator fg of the CAM model based on the LLNL observations. Column 1: distance 
metric d used to define the Lipschitz constant. Columns 2~4'- Binomial lower confidence bounds for quartiles 
of the pointwise error, obtained by inverting binomial tests. Column 5: 95% lower confidence bound for the 
integral of the pointwise error over the entire domain [0, l] p , based on inverting a z-test. Columns 2-5 are 
expressed as a fraction of K /2 — half the empirical Lipschitz constant. Results are based on 10,000 uniform 
random samples from [0, l] p . 



5. Extensions. The theory developed in this paper focuses on the maximum uncertainty 
over all w G [0, l] p . We have argued that this maximum uncertainty is important in some 
applications. In other applications, whether the uncertainty is anywhere large might be less 
interesting than the fraction of the domain [0, l] p on which the uncertainty is large. 

It is possible to estimate the fraction of [0, l] p for which e* exceeds any fixed threshold 
e > by sampling. By drawing values w 6 [0, l] p at random and evaluating e* at each selected 
point, it is possible to construct a lower confidence bound for the fraction of the domain [0, l] p 
for which the potential error exceeds any given threshold, and to construct confidence bounds 
for quantiles of the potential error. 

Table 5.1 shows the results for the CAM simulations, based on 10,000 random samples 
from [0, l] p . Even the lower quartiles are a large fraction of K. For instance, at confidence 
level 95%, the potential error under the sup-norm metric exceeds 71.7% of K/2 on at least 
50% of the domain. 

6. Conclusions. We find a lower bound on the minimum (over emulators) maximum (over 
regular functions that agree with the observations) error of emulators of a function / based 
on n observations. This "mini-minimax" uncertainty is optimistic because it assumes that / 
is no less regular (i.e., not "rougher") than the n observations require / to be. It measures 
the potential error of the best emulator of / at the point w in the domain of /: there are 
functions g and h that are no less regular than /, that agree with / at the n observations, and 
that differ at the point w by twice this potential error. No emulator can approximate both 
functions with less error than the potential error at w. The observations and the regularity 
condition do not rule out g or h: the true / could be either. 

In some problems, every emulator based on any tractable number of observations of / has 
large maximum potential error (and the potential error is large over much of the domain), 
even if / is no less regular than it is observed to be. That is, there are functions g and h that 
agree perfectly with the observations, are no less regular than the observations require them 
to be, and yet differ by a large amount at some point in the domain of /. When that occurs, 
all emulators of / are potentially substantially incorrect. 

We give sufficient conditions under which all emulators are potentially substantially in- 
correct. The conditions depend only on the observed values of /; they can be computed from 
the same observations used to train an emulator, at a cost that typically is small compared 
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with the cost of generating those observations. The conditions are sufficient but not necessary, 
because / could be less regular than any finite set of observations reveals it to be. It is not 
possible to give necessary conditions that depend only on the observed values of / (a priori 
bounds on the regularity of / would be needed). 

Our conditions seem likely to hold for many HEB applications of societal interest, such 
as climate modeling: the curse of (-High) dimensionality makes it necessary to observe / at 
many points, unless it is known a priori that / is extremely regular; the .Expense makes it 
intractable to observe / at that many points; and the i?lack-box nature of / means we do not 
typically know a priori how regular / is. Indeed, we show quantitatively that these conditions 
hold for a common measure of regularity and a large climate modeling dataset. 

When the maximum potential error in approximating / everywhere by a constant — the 
value of / at the center of the domain — is no larger than the maximum potential error in 
approximating / from any tractable number of observations, emulators may not be useful. In 
that circumstance, no emulator can reliably model / as a function of its input w € [0, l] p . 

Common techniques for validating emulators do not reveal whether emulators are poten- 
tially substantially incorrect, either because they make strong assumptions about / that are 
based neither on the observations nor on real knowledge of the properties of / (e.g., Bayesian 
approaches), or because they estimate the average error with respect to some measure, rather 
than the potential error (e.g., validation using held-out observations). Our results reflect 
the long-standing tension between worst-case analysis and average-case analysis, 6 between 
Bayesian and minimax measures of uncertainty, and between simple, high-bias models and 
complex models that are harder to estimate and interpret. The complexity of many HEB mod- 
els may be on the far side of the tradeoff, necessitating the use of emulators that are destined 
to be highly uncertain because computational constraints preclude adequate sampling. 

Reducing the potential error of emulators in HEB problems requires either more informa- 
tion about / (knowledge, not merely assumptions 7 ), or changing the measure of uncertainty — 
changing the scientific question. Both tactics are application-specific: the underlying science 
dictates the conditions that actually hold for / and the senses in which it is useful to approx- 
imate /. 

Finally, it is not clear that emulators help address the most interesting questions. Ap- 
proximating / pointwise is not typically the ultimate goal; in many problems, most prop- 
erties of / are nuisance parameters. It could well be that the important questions about / 
can be answered more directly, without estimating / as a whole. For example, for global 
optimisation — finding maxima or minima — a form of adaptive sampling known as multi-start 
methods yield good results [21]. It may be possible to develop an adaptive approach to finding 
level sets — another common use for emulators — too [22]. 

There may not always be a remedy, however. Some research questions simply cannot be 
answered through simulation at present. Employing complex emulators and massive compu- 

6 However, as section 5 shows, the average potential error for the CAM model is also quite large. 

'Common additional conditions include the following: parameters have only low order interactions; the 
second derivative has an upper bound; the third derivative has a limited number of knots; the integral of the 
squared derivative of the model is bounded [20]. There are problems in which conditions like these may reflect 
actual knowledge about /. However, such conditions tend to be difficult to verify: simulation is perhaps most 
valuable when the underlying equations are not amenable to mathematical analysis. 
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tational resources can distract us from this reality. With greater scientific understanding of 
the underlying function, it may become possible to reduce uncertainty in HEB simulations 
to useful levels. Until then, simpler models — with caveats about what they omit — may be a 
better basis for scientific inference. 
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Appendix A. Proofs. 

To start, we relate the maximum potential error to an intersection of intervals in the range 
of /. For real x an d p, define the interval (one-dimensional ball) 



i(x,p) 



[x- p,x + p], P > o 

0, otherwise. 



If / is an interval, p(I) denotes its length; for instance, p(I(x,p)) = max(0, 2p). 
For g eC[0, 1} P and kGK, 



1 



2 

- <j min {g(x) + Kd(x, w)} — max {g(x) — nd(x, w)} 

2 I xdX x&X 



(A.l) 



2^ 



max {g{x) — nd(x, w)} , min {g{x) + nd{x, w)} 



\xex / 

Lemma 7. Fix a € [0, 1], pi, . . . , p n G [0, oo) and xi, ■ ■ ■ , Xn GK. Let I\ = HILi HXi> Pi) an d 
I a = HXi, api). Then ap (A) > /U (J a ). 

Proof. Because the intersection of intervals is itself an interval, there exist xo an d po satisfying 

I a = I(xo,po)- 



Fix i G 1, . . . , n. Then 
It follows that 
Then 

Because a < 1 and pi > 0, 
Finally, 



I(xo,Po) C I(xi,api 



Xo- Po>Xi~ api. 



a \Pi - —) >Xi~ Xo- 
>Xi~ Xo- 



a 
Po 



a 



Po . 

XO >Xi~ Pi- 

a 



By symmetric reasoning we also have that 



Therefore, 



. Po - 

Xo H < Xi + Pi- 

a 



l(xoA)cI( X i,Pi). 
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Because i was arbitrary, 
Hence, 



H(Ii)>H[I [Xo,— = " 



□ 

Lemma 7 is used in the proof of Theorem 3, below. 
Proposition 1. supe* = £*. 

Proof. 

Step 1: e+ and e~ are continuous. 

For fixed x G X, if g £ C[0, l] p , g(x) ± K,d(x,w) is continuous in u;. Because the minimum or 
maximum of a finite set of continuous functions is also continuous, mm x£ x [9( x ) + Kd(x, w)] 
and max l6 x [g(x) — nd(x, w)} are continuous in w too. Thus e+ and e~ are continuous func- 
tions on [0, l] p . 

Step 2: e+ and e~ agree with f on X. 

Suppose g has Lipschitz constant k and g agrees with / on X. Then for any x € X and 
w€[0,lp, 

f(x) — K,d(w, x) < g{w) < f(x) + nd(w, x). 



Hence, 



Ifw&X, 



max [f(x) — nd(x, w)] < g(w) < min [f(x) + nd(x, w)] . 



min [f(x) + nd(x, w)] = max [/(x) — nd(x, w)] = f(w). 



Thus, e+ and e K € T K agree with / on X. 

Step 3: e+ and e~ have Lipschitz constant k. 
For v, w € [0, l] p , 3x, y £ X satisfying 

e t( v ) = f( x ) + Kd(x,v) and e+(w) = f(y) + nd(y,w). 

Suppose without loss of generality that e + (v) > e + (w). By construction, e^(v) < f{y) + 
nd(y,v). Hence 

< e+(«) - e+O) < f{y) + nd(y, v) - e+O) 

= f{y) + Kd{y, v) - f(y) - Kd(y, w) 
< n(d(y,v) - d(y,w)) 

<Kd(v,y), (A.2) 

by the triangle inequality. Hence e+ has Lipschitz constant k. An analogous argument shows 
that e~ also has Lipschitz constant k. Combining this with steps 1 and 2 shows that e + and 
e~ are in T K . 
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Step 4-' e K is the pointwise minimum of T K and e+ is the pointwise maximum of T K . 
Suppose to the contrary that there exists w 6 [0, l] p , and g G T K for which 

g(w) > f(x) + nd(x, w). 

Recall that g G T K implies that g(x) = f{x) Vx G X. Hence 

g{w) — g{x) > f{x) + nd(x, w) — f(x) = nd(x, w). 

That is, g has a Lipschitz constant greater than k, a contradiction. Hence, e~^{w) = sup{g(w) : 
g G J- K } for all w G [0, The same argument, mutatis mutandi, shows that e~(w) = 
inf{3(-u;) : g G F K } for all w G [0, If. 

Step 5: sup e* = £* . 
Now 

sup e*(» = sup - [e+(ifl)-e~f»] 

u>G[0,1]p «;G[0,1]p z 

= - sup < sup g(w) — inf h{w) > 

= \ sup i sup - h{w)\ \ 

A u>6[o,i] p [gAe^K J 

= ^ sup sup \g(w) - h(w)\ 
1 g,h&T K i/jg[0,1]p 

= ^sup{||5f - /j||oo : g,h G 



□ 



Theorem 3. For any X, if supe*-, > A A, then £x(f) > A A. 
Proof. Let = arg max,,, e*^ (w) . Then 

= supe^ (A.3) 

> e* K (w*) 

>f-4K) ( A - 4 ) 

> 5- • AA (A.5) 
"if 

= A A'. 

(A.3) is a consequence of proposition 1. (A.5) follows from (A. 4) by hypothesis. (A. 4) is 

a consequence of lemma 7: Let a = K/K < 1. For, i = 1, . . . ,#X, let = /(xj) and 

Xi = Kd(x,w). Then, by (A.l), fi(h)/2 = e K and n{I a )/2 = eg- □ 
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Lemma 5. Let = (0, . . . , 0) and 1 = (1, . . . , 1) . For v, w € [0, l] p , 

d(v, w) < d(v) = max (d(v, 0),d(v, 1)) . 
Proof. Let wr^ denote the i th component of w. Then 



d(v, w) 



max \V(i) — W(j\ 



< max max \v 
ie{i,-,p} <5e{o,i} 



(0 



max max va\ — va\ 
ie{i,..., P } y e{o,l} 1 W 

max max I vm — I 

max d(t>, y) 
2/e{0,i} 

max(d(v, 0), 1)). 



Proposition 6. 



1 

supe^ < - 



Proof. Fix w £ [0, If. Then, 
* r \ 1 



mm 



f(x) + £d(x) 



max 



f( x ) - kd(x) 



fl I[f(x),Kd( 

= - < min \f(x) + Kd(a 
2 I x£X L 



max 



/(x) - Kd(x) 



□ 



(A.6) 
(A.7) 



where (A.6) follows from (A.l), and (A.7) follows from lemma 5. Because the right-hand side 
of this inequality does not depend on w, the proposition follows by taking suprema. □ 



